﻿using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;

namespace RGeos.Framework.OTL
{
    /// <summary>
    /// 导线荷载和张力计算
    /// </summary>
    public class OTLStrainCaculate
    {
        //重力加速度
        const double g = 9.80665;
        public static double JFC(Double a, Double b)
        {
            double a1, t, T1;
            int C;
            double vJFC;
            C = Math.Sign(a);
            if (C != 0)
            {
                a1 = Math.Abs(a);
                t = 13.5 * b / Math.Pow(a1, 3) - C;
                if (t > 1)
                {
                    T1 = Math.Log(t + Math.Sqrt(t * t - 1));
                    vJFC = (2 * Math.Cosh(T1 / 3) - C) * a1 / 3;
                }
                else
                {
                    T1 = Math.Acos(t);
                    vJFC = (2 * Math.Cos(T1 / 3) - C) * a1 / 3;
                }
            }
            else
            {
                vJFC = Math.Pow(b, (1 / 3));
            }
            return vJFC;
        }
        /// <summary>
        /// 计算荷载
        /// </summary>
        /// <param name="V">风速</param>
        /// <param name="b">冰厚</param>
        /// <param name="W">单位长度重量</param>
        /// <param name="d">直径</param>
        /// <returns></returns>
        public static double GetP(double V, double b, double W, double d)
        {
            double g3 = g * W + 0.0009 * g * Math.PI * b * (b + d);
            double g5 = g * Getk(d, b) * GetC(V) * (d + 2 * b) * Math.Pow(V, 2) / 16000;
            return Math.Sqrt(Math.Pow(g3, 2) + Math.Pow(g5, 2));
        }
        /// <summary>
        /// 计算风荷载
        /// </summary>
        /// <param name="V">风速</param>
        /// <param name="b">冰厚</param>
        /// <param name="w">导线单位长度重量</param>
        /// <param name="d">导线直径</param>
        /// <returns></returns>
        public static double GetP5(double V, double b, double w, double d)
        {
            double p5 = g * Getk(d, b) * GetC(V) * (d + 2 * b) * Math.Pow(V, 2) / 16000;
            return p5;
        }
        /// <summary>
        /// 风载体型系数
        /// </summary>
        /// <param name="d">导线直径</param>
        /// <param name="b">冰厚</param>
        /// <returns></returns>
        static double Getk(double d, double b)
        {
            double Getk;
            if (b != 0)
            {
                Getk = 1.2;
            }
            else if (d >= 17)
            {
                Getk = 1.1;
            }
            else
            {
                Getk = 1.2;
            }
            return Getk;
        }
        /// <summary>
        /// 风速不均匀系数
        /// </summary>
        /// <param name="V">风速</param>
        /// <returns></returns>
        static double GetC(double V)
        {
            if ((V >= 30) && (V < 35))
            {
                return 0.75;
            }
            if (V >= 20 && V < 30)
            {
                return 0.85;
            }
            if (V < 20)
            {
                return 1;
            }
            if (V >= 35)
            {
                return 0.7;
            }
            return 0;
        }
        /// <summary>
        /// 计算临界档距
        /// </summary>
        /// <param name="EA">E*A--弹性系数乘以截面积</param>
        /// <param name="EA1"></param>
        /// <param name="r1">四种工况数组</param>
        /// <returns>临界档距，第一列为临界档距值，第二列为控制的气象条件在数组中的索引</returns>
        public static double[,] LJDJ(double EA, double EA1, double[,] r1)
        {
            double[] F, P, t;
            F = new double[4];
            P = new double[4];
            t = new double[4];
            int i, j, mini, maxj;
            double[] Lcr = new double[7];
            double[] Lcr1 = new double[7];
            double[,] kz = new double[4, 2];
            double[,] result = new double[4, 2];

            //EA = Application.Worksheets("sheet1").Range("O2").Value
            //'EA1 = Application.Worksheets("sheet1").Range("P2").Value
            //'Set r1 = Application.Worksheets("sheet1").Range("E8:G11")
            for (i = 0; i < 4; i++)
            {
                t[i] = r1[i, 0];
                P[i] = r1[i, 1];
                F[i] = r1[i, 2];
            }
            //'Set r1 = Nothing
            int n = 0;
            for (i = 0; i < 3; i++)
            {
                for (j = i + 1; j < 4; j++)
                {
                    double l = ((double)24 * (F[i] - F[j]) + 24 * EA1 * (t[i] - t[j])) / EA;
                    //l = l / (((P[i] / F[i]) * (P[i] / F[i]) - P[j] / F[j]) * ((P[i] / F[i]) * (P[i] / F[i]) - P[j] / F[j]));
                    l = l / (Math.Pow((P[i] / F[i]), 2) - Math.Pow((P[j] / F[j]), 2));
                    if (l > 0)
                    {
                        Lcr[n] = Math.Sqrt(l);
                        n = n + 1;
                    }
                }
            }
            Lcr[n] = 0;
            //n = n - 1;

            for (i = 0; i < n - 1; i++)
            {
                mini = i;
                for (j = i + 1; j < n; j++)
                {
                    if (Lcr[j] < Lcr[mini])
                    {
                        mini = j;
                    }
                }
                double tmp = Lcr[i];
                Lcr[i] = Lcr[mini];
                Lcr[mini] = tmp;
            }

            Lcr1[0] = Lcr[0] / 2;
            for (i = 1; i < n; i++)
            {
                Lcr1[i] = (Lcr[i] + Lcr[i - 1]) / 2;
            }
            if (n < 1)
            {
                return result;
            }
            Lcr1[n] = Lcr[n - 1] + 10;

            int n1 = 0;
            for (i = 0; i < n + 1; i++)
            {
                maxj = 0;
                for (j = 1; j < 4; j++)
                {
                    double FX1 = EA * Math.Pow((P[j] * Lcr1[i] / F[j]), 2) / 24 - (F[j] + EA1 * t[j]);
                    double FX2 = EA * Math.Pow((P[maxj] * Lcr1[i] / F[maxj]), 2) / 24 - (F[maxj] + EA1 * t[maxj]);
                    if (FX1 > FX2)
                    {
                        maxj = j;
                    }
                }
                if (n1 == 0)
                {
                    kz[n1, 0] = Lcr[i];
                    kz[n1, 1] = maxj;
                    n1 = n1 + 1;
                }
                else
                {
                    if (kz[n1 - 1, 1] == maxj)
                    {
                        kz[n1 - 1, 0] = Lcr[i];
                    }
                    else
                    {
                        kz[n1, 0] = Lcr[i];
                        kz[n1, 1] = maxj;
                        n1 = n1 + 1;
                    }
                }
            }

            for (i = 0; i < 4; i++)
            {
                result[i, 0] = kz[i, 0];
                result[i, 1] = Math.Ceiling(kz[i, 1]);
            }
            return result;
        }
        /// <summary>
        /// 计算张力临界值
        /// </summary>
        /// <param name="l">代表档距</param>
        /// <param name="tc">温度</param>
        /// <param name="Pc">综合载荷</param>
        /// <param name="r1">临界档距，5*3数组</param>
        /// <param name="r2">四种工况数组</param>
        /// <returns>水平张力，单位为N</returns>
        public static double GetT(double l, double tc, double Pc, double[,] r1, double[,] r2)
        {
            double EA = r1[0, 0];
            double EA1 = r1[0, 1];
            int km = 0;
            int i = 0;
            for (i = 1; i < 4; i++)
            {
                if (r1[i, 0] == 0)
                {
                    km = i - 1;
                    break;
                }
            }
            int k = 0;

            for (i = 1; i <= km; i++)
            {
                if (l <= r1[i, 0])
                {
                    k = (int)r1[i, 1];
                    break;
                }
            }
            if (k == 0)
            {
                k = (int)r1[i, 1];
            }

            double Tk = r2[k, 0];
            double Pk = r2[k, 1];
            double Fk = r2[k, 2];


            double x = (Pk * l) * (Pk * l) * EA / (24 * Math.Pow(Fk, 2)) - Fk + EA1 * (tc - Tk);
            double y = (Pc * l) * (Pc * l) * EA / 24;
            return JFC(x, y);
        }

        public static double GetPolylineK(double ZHZH, double A, double zL)
        {
            double K0 = ZHZH / zL / 8;
            return K0;
        }
    }
}
